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Abstract 

Dynamics of a coarse-grained model for the room-temperature ionic liquid, l-ethyl-3- 
methylimidazolium hexafluorophosphate, couched in the united-atom site representation are 
studied via molecular dynamics simulations. The dynamically heterogeneous behavior of the 
model resembles that of fragile supercooled liquids. At or close to room temperature, the 
model ionic liquid exhibits slow dynamics, characterized by nonexponential structural relax- 
ation and subdiffusive behavior. The structural relaxation time, closely related to the viscosity, 
shows a super-Arrhenius behavior. Local excitations, defined as displacement of an ion ex- 
ceeding a threshold distance, are found to be mainly responsible for structural relaxation in the 
alternating structure of cations and anions. As the temperature is lowered, excitations become 
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progressively more correlated. This results in the decoupling of exchange and persistence 
times, reflecting a violation of the Stokes-Einstein relation. 

1 Introduction 

Due to their important potential as environmentally benign alternatives to conventional toxic or- 
ganic solvents, room-temperature ionic liquids (RTILs) have attracted considerable attention re- 
cently.- - - According to theoretical- - — and experimental— - — studies on solvation and rotational 
dynamics of RTILs, their long-time behaviors are characterized by nonexponential decay. This 
implies that RTILs are dynamically inhomogeneous and their local relaxation is widely distributed 
in time and space. 510 The clustered mobile and immobile ions observed in recent molecular dy- 
namics (MD) simulation studies are ascribed to inhomogeneous dynamics in RTILs. — ^ 

Dynamic heterogeneity often invoked to explain many peculiar properties of supercooled liq- 
uids 19,20 refers to the enhanced temporal correlation of their local dynamic states with a decrease in 
temperature. From the viewpoint of facilitated motions, dynamics of supercooled liquids are domi- 
nated by fluctuations.— & According to several studies based on lattice models, called kinetically- 
constrained models (KCMs),— ^ non-trivial structures in the space-time trajectory arising from 
dynamic constraints in the KCM description accurately reproduces many of the dynamical prop- 
erties of supercooled liquids. 2 ^ 5 - At the molecular level, it is found that trajectories of individual 
particles in atomistic models of supercooled liquids are in general governed by dynamic fluctua- 
tions and thus cannot be predicted from static properties, such as structures. 26 Meanwhile, recent 
studies have attempted to correlate length-scale dependent heterogeneous dynamics with liquid 
structures on the basis of the dynamic propensity calculated from the isoconfigurational ensem- 
ble. 27 Despite these efforts, the origin of persisting dynamic correlations and the potential link of 
the dynamic correlations to structures still remain open questions in understanding dynamics of 
glassy liquids. 

Glassy dynamics of supercooled liquids are characterized by many unique features such as the 
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stretched exponential decay of time correlation functions, subdiffusive behavior in the intermediate 
time scale, decoupling of exchange and persistence times, 22 and breakdown of the Stokes-Einstein 
(SE) relation.— A variety of models of supercooled liquids, e.g., binary Lennard- Jones (LJ) mix- 
ture, 29 supercooled water, 30 Weeks-Chandler- Anderson (WCA) mixture, 31 and KCMs, 2 -' 3 - 2 -' 3 - 3 - re- 
produce these interesting features, which are generally believed to be intimately related to the 
dynamic heterogeneity. This seemingly universal nature of dynamic heterogeneity and related 
phenomena is a main motivation to study glassy dynamics of ionic liquids, the molecular details 
and interaction potentials of which are quite different from those of supercooled liquids. Specif- 
ically, strong electrostatic interactions in RTILs lead to an alternating structure of cations and 
anions, which is believed to generate a significant memory effect on their dynamics. 10 Since the 
cage structure would exert a substantial influence on correlations between local excitation events, 
it would be both interesting and important to see if dynamic behaviors of RTILs bear any resem- 
blance to those of supercooled liquids. 

One of the main challenges in studying dynamics of RTILs (and more generally glassy liquid 
systems) is the time scale. As is well known, many glassy systems exhibit extremely slow dynam- 
ics, in particular, at low temperatures. Therefore it is very difficult, if not impossible, to probe 
their long-time behaviors directly via atomistic MD simulations. However, it is precisely these 
dynamics at long times and their variations with temperatures that are of especial interest. We thus 
employ a coarse-grained model to make simulations efficient and analyze long-time dynamics. 
We note that this approach is not new. Previous efforts 34-36 based on similar descriptions have 
provided useful insights into RTIL properties, both structure and dynamics. In the present paper, 
we focus on the characterization of glassy dynamics and dynamic heterogeneity of RTILs. As a 
prototypical RTIL, we study l-ethyl-3-methylimidazoliumhexafluorophosphate (EMI + PF6 ). 

To understand translational dynamics of EMI + PF 6 , we analyze the self-intermediate scatter- 
ing function and mean square ion displacements of our model RTIL system. We also examine the 
temperature dependence of its structural relaxation time and related fragility. 37 Our results sug- 
gest that EMI + PF 6 is a fragile glass former and it violates the SE relation at low temperatures. 
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We regard correlations between local dynamic states as an essential aspect of the aforementioned 
dynamic heterogeneity of glassy liquids. Accordingly, we investigate distributions of the persis- 
tence and exchange times, i.e., waiting times for the first and subsequent excitations, by monitoring 
the motions of individual ions. To understand the influence of Coulombic interactions on glassy 
dynamics, we also make contact with non-ionic model systems of supercooled liquids. 

The outline of this paper is as follows: In Sec.[2l we introduce the coarse-grained model of an 
ionic liquid and briefly describe the simulation methods. In Sec. [3l we present the MD results 
on glassy dynamics, including fragility, nonexponential relaxation, subdiffusion, and violation 
of the SE relation. The correlation of local events together with the decoupling of persistence 
and exchange times is demonstrated in Sec. HI while the roles of the Coulombic interactions are 
discussed in Sec. [51 Concluding remarks are offered in Sec. [6] 

2 RTIL models 

In this section, we give a brief explanation of our coarse-grained model 34 for EMI + PFg and com- 
pare its results for structure and translational dynamics with those of a more atomistic description. 

2.1 Coarse-grained model 

Our coarse-grained model is based on the RTIL description used by Kim and coworkers— to study 
various processes in RTILs.- 1 ^^ Specifically, they employed the united atom representation for 
the methyl group (Ml) as well as the CH2 (El) and CH3 (M3) moieties of the ethyl group of 
cations. They used the AMBER force field 40 for the Lennard- Jones (LJ) interactions and partial 
charge assignments 41 proposed by Lynden-Bell and coworkers for Coulombic interactions. PF6 
was also described as a united atom. Hereafter, this model will be referred to as the AM description. 
In our coarse-grained model (CGM), we further simplify the cation description by representing the 
imidazole ring and H atoms directly attached to it as a single atomic site Tl (Fig. [Q). Thus each 
EMI + ion consists of 4 united-atom sites, Ml, M3, El and Tl in CGM. The LJ parameters of 
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Figure 1: Coarse-graining scheme. EMI + is reduced to a 4-atom cation, where 5 atoms of the 
imidazolium ring of EMI + and 3 hydrogen atoms directly attached to it are collapsed to a single 
united atom (Tl) in the coarse-grained model. The methyl and ethyl groups of EMI + are repre- 
sented, respectively, as one-atom (Ml) and two-atom (M3 and El) moieties and PF 6 is simplified 
as a united atom as in Ref. 38. 
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Table 1 : The LJ parameters, partial charges, and masses of coarse-grained atoms 



atom 


on (A) 


£„ (kJ/mol) 


qi (c) 


mass (amu) 


Ml 


3.905 


0.7330 


0.316 


15.04092 


Tl 


4.800 


1.5000 


0.368 


67.08860 


El 


3.800 


0.4943 


0.240 


14.03298 


M3 


3.800 


0.7540 


0.076 


15.04092 


PF6 


5.600 


1.6680 


-1.000 


144.97440 



Tl were adjusted, so that CGM reproduces the liquid structure of the atomistic AM description 
reasonably well (see below). For Ml, M3, El and PF 6 , we used the parameters of AM without 
any further adjustment. 

We performed MD simulations of EMI + PF 6 in both the AM and CGM representations using 
the DL_POLY program. 42 Atoms i and j at positions r, and Tj interact with each other through the 
LJ and Coulomb potentials: 



U tj = 4 £i 



'ij 



12 



'tj 



. HiH] 



(1) 



where r,- y = |r, — r ; | is the distance between the two atoms. The parameters of CGM employed in 
the present study are compiled in Table CD For the AM parameters, the reader is referred to Ref. 
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. For LJ interactions between unlike atoms, the Lorentz-Berthelot combining rules were used. 43 
The simulation cell of the CGM ionic liquid comprises 512 pairs of rigid cations and an- 
ions. We performed simulations in the canonical ensemble at six different temperatures, T = 
300, 350, 400, 475, 600 and 800 K, using the Nose-Hoover thermostat and at density p = 1 .3 1 g/cm 3 . 
Periodic and cubic boundary conditions were employed and long-range electrostatic interactions 
were computed via the Ewald method. Starting from a crystal configuration obtained by alternat- 
ing the cations and anions, we equilibrated the system for 2 ns prior to production runs at 800 K. 
At lower temperatures, we used as an initial configuration one of the equilibrated configurations at 
higher T that is closest to the temperature of the system under consideration and equilibrated the 
system for 2 ns at 600 K and 475 K, 5 ns at 400 K, and 10ns at 350 K and 300 K. Production runs 
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r (A) r (A) 

Figure 2: Radial distributions of ions in EMI+PF 6 at (a) 350 K and (b) 400 K. The results of CGM 
and AM are plotted in solid and dashed lines, respectively. The center of mass position is used to 
represent cation locations in CGM and AM. 

following equilibration were 5 ns in length for 800 K and 600 K, 10 ns for 475 K and 400 K, 50 ns 
for 350 K, and 60 ns for 300 K. At each thermodynamic condition, we carried out six independent 
production runs, from which the averages were computed. Thus we used, for example, six 60ns 
trajectories to analyze various dynamic quantities at 300 K in CGM. 

In the case of AM, we considered 112 pairs of EMI + and PFg ions with p = 1.31 g/cm 3 at 
T = 350 and 400 K. At either temperature, we performed three independent simulations, each of 
which was carried out with 10ns equilibration, followed by a 40ns trajectory. Ensemble averages 
were calculated using three trajectories thus obtained. 

2.2 Structure 

Here we consider MD results of the CGM and AM for structure to gain insight into how realistic the 
former description is. In Fig.[2l their results for radial distributions of ions at T = 350 K and 400 K 
are compared. There we employed the center of mass to represent the cation positions for CGM 
and AM. We notice that CGM captures the RTIL structure of AM very well. The two models yield 
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Figure 3: Mean square displacement A(t) of (a) cations and (b) anions in CGM and AM at T = 
350 K and 400 K. 

a excellent agreement in both the peak positions and heights, including minor secondary peaks, 
e.g., structure near 7 A in the cation-anion distribution. Even in the case of the main peak of the 
cation-anion distribution which shows the largest deviation between the two, the discrepancy in 
peak location is only ~ 0.02 A. Considering the drastic approximation of the planar imidazole ring 
as a spherical atom, we think that overall the coarse-grained model fares very well in reproducing 
the RTIL structure of AM. 



2.3 Translational dynamics 

To understand the effect of our coarse-graining (Fig.Q]) on system dynamics, we examine ion trans- 
lational motions by considering their mean square displacement, A(f) = (N^ l Y. ! iLi \ r i(t) ~ r /'(0) | 2 )' 
and self-intermediate scattering function, F s (qo,t) = (N^ 1 T!jLi e ! l '[ r 'W _r '(°)]}. Here (■ • • ) denotes 
the equilibrium ensemble average, N the number of ions, and qo the wave vector corresponding to 
the position of the first peak in the static structure factor. 

MD results in Figs. [3] and |4] show that the coarse-grained model exhibits subdiffusive behavior 
and nonexponential relaxation, consistent with the atomistic model. We, however, notice that 
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Figure 4: Self-intermediate scattering function F s (qQ,t) for (a) cations and (b) anions in CGM and 
AM at T = 350 K and 400 K. The wave vector employed in (a) and (b) are, respectively, go = 
0.858 A and 0.878 A , corresponding to the positions of the first peak in their static structure 
factor. 

Table 2: MD results for the diffusion constant D and the structural relaxation time x a for CGM 



r(K) 


D(A"/ps) 






(ps) 




cation 


anion 




cation 


anion 


300 


5.4x 10~ 5 


9.6 x 10 


-6 


59600 


244000 


350 


4.3 x 10~ 4 


9.7 x 10 


-5 


5030 


12700 


400 


2.2 x 10~ 3 


7.2 x 10 


-4 


779 


1610 


475 


8.8x 10~ 3 


4.0 x 10 


-3 


154 


261 


600 


3.3 x 10~ 2 


2.0 x 10- 


-2 


36.1 


50.2 


800 


9.9x 10~ 2 


6.8 x 10- 


-2 


11.4 


13.5 



dynamics of the former proceed faster than the latter. To quantify this, we calculated the diffusion 
constant D = lim f ^ 00 [6(? — ?o)] _1 (W -1 Y.?=i [ r /(0 — r K ? o)] 2 ). where to denotes the time when the 
Fickian behaviors appear in A(t), and determined the structural relaxation time x a via 

F s (qo,T a )=e- 1 . (2) 
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The MD results for D and x a are compiled in Tabled At 400 K, we obtained Z) = 2.2 x 10 3 A /ps 

° 2 

and 7.2 x 10~ 4 A /ps for cations and anions, respectively, for CGM, while the corresponding val- 

a 2 9 2 

ues for AM were 6.5 x 10" 4 A /ps and 1.4 xlO~ 4 A /ps. As for T a , the coarse-grained model 
yields 779 ps and 1610ps for cations and anions, respectively, whereas the AM description results 
in significantly longer relaxation times, 3680 ps and 9010 ps. Because of this discrepancy, i.e., 
CGM is faster than AM in dynamics, A(t) and F s (q,t) of AM at 400 K are in better accord with 
those of CGM at 350 K than at 400 K (Figs. [3]and|4]). (For CGM at 350 K, simulations yield 

9 2 c ° 2 

D = 4.3 x 10~ 4 A /ps and 9.7 x 10~ 5 A /ps, and x a = 5030ps and 12700ps for cations and anions, 
respectively.) This is not surprising in that structure and dynamics of liquid systems, including 
ionic liquids, depend on molecular shape (i.e., LJ interactions) and charge distributions (viz., elec- 
trostatic interaction) of constituent particles.— 1 ^ In addition to the difference in cation geometry, 
negative partial charges of nitrogen atoms of EMI + present in the AM description are completely 
absent in CGM because of the united atom representation of the imidazole ring. This simplifi- 
cation reduces charge anisotropy of cations and thus frustration in the structure and dynamics of 
CGM. We therefore expect that both rotational and translational dynamics would be accelerated in 
CGM, compared to AM. Here we consider only the translational dynamics of cations and anions 
and postpone the rotational dynamics for a future study. 



3 Glassy dynamics 

In this section, we analyze the glassy behavior of the coarse-grained RTIL with the aid of MD 
simulation results. We examine structural relaxation and ion diffusion at various temperatures and 
demonstrate that our model belongs to the class of fragile glass formers and violates the Stokes- 
Einstein relation. 
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Figure 5: Self-intermediate scattering function F s (qo,t) for (a) cations and (b) anions in a coarse- 
grained ionic liquid at various temperatures. The wave vector qo is the same as in Fig. @] In 
the a -relaxation regime, F s (qo,t) exhibits nonexponential decay, which is well described by 



3.1 Structural relaxation and fragility 

The structural relaxation through ion translational dynamics is usually described by the self- 
intermediate scattering function F s (qo,t). In Fig. |51 the CGM predictions for F s (qo,t) at six dif- 
ferent temperatures are displayed. At high T, the ionic liquid behaves almost like a normal liquid; 
the thermal fluctuations dominate over the constraints of caged structures of ions, mainly formed 
by counterions. As the temperature decreases, characteristics of glassy dynamics, such as subdif- 
fusivity and nonexponential relaxation, become more apparent. 

At low T, the presence of the plateau regime (/3 relaxation) and slow a relaxation, which are 
hallmarks of supercooled liquids, is rather prominent. The contribution of inertial dynamics in 
the first 0.1 ps or so accounts for less than 10% of the entire relaxation of F s (qo,t), and structural 
correlation persists for several decades in time, thereby indicating the highly viscous RTIL envi- 
ronment. The slow nonexponential relaxation subsequent to the plateau regime is well described 
by a stretched exponential function, cexp[— (?/tq)^]. At 300 K, the exponent /3 is found to be 0.64 



c exp [—(t/TQ)P] for all temperatures . 
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Figure 6: Temperature dependence of structural relaxation time x a , Eq.|2l Dashed and dotted lines 
are the fits to MD results using (a) x a = c\ exp(d\/T 2 ) and (b) x a = C2exp(d2/T). Error bars 
represent the maximum and the minimum values of % a among six independent trajectories. 

and 0.59 for cations and anions, respectively. A substantial deviation of these /3 values from unity 
is another good indicator of the glassy dynamics in the RTIL system. As T increases, so does /3. 
At 800 K, the highest temperature we studied, the /3 values are 0.89 and 0.92. Even though greatly 
enhanced thermal fluctuations at this temperature accelerate structural relaxation immensely by 
more than four orders of magnitude compared to that of 300 K, F s (qo,t) still maintains a nonexpo- 
nential character, presumably due to its high pressure condition. We thus expect that if we lower 
the pressure by reducing its density, the structural relaxation would become a single-exponential 
decay. 

We turn to the temperature dependence of % a in the CGM description presented in Figure [6] 
and Table [2l Since slow structural relaxation at 300 K does not allow the determination of its 
% a directly from the simulation results, we estimated it by employing a stretched exponential fit 
to F s (qo,t). At all other temperatures, x a was determined using the MD results for F s (qo,t) in 
Eq. [2l The most salient aspect of our results in Figure [6] is that x a does not follow the Arrhenius 
law % a 06 exp(<i2/7 1 )- Rather the structural relaxation time shows a super- Arrhenius behavior; 
specifically, it varies with the temperature as x a exp(Ji/r 2 ). This means that the CGM RTIL 
studied here resembles a fragile glass former in the temperature dependence. 
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For clarity, we make a couple of remarks here. First, while the RTIL density is assumed 
to be fixed in the present study, it tends to decrease with increasing temperature for real ionic 
liquids. This density variation, if incorporated, would lead to acceleration of structural relaxation 
in CGM at high T, compared to the results in Fig. [6l This would in turn strengthen the super- 
Arrhenius character of x a and thus enhance the fragile behavior of CGM. Second, as pointed out 
in Sec. 12.31 above, temperatures of CGM and AM do not coincide. In other words, the temperature 
of CGM does not correspond to the actual temperature of the atomistic system. As a consequence, 
our finding that EMI + PF 6 is a fragile glass former based on the CGM description might not be 
directly applicable to the real ionic liquid. We however note that a super- Arrhenius temperature 
dependence was observed in recent measurements in a similar RTIL. 44 We thus believe that our 
result on the fragility of EMI + PFg is robust. 

3.2 Breakdown of the Stokes-Einstein Relation 

In a normal liquid, the diffusion constant D is usually related to the viscosity r\ via the Stokes- 
Einstein relation 

T 

Doc - . (3) 
V 

For convenience, we consider another relation 

Doc — , (4) 

which is equivalent to Eq. [3] if the structural relaxation time is proportional to r}/T. Eq. 0] results 
when translational motions of constituent particles are described by a normal diffusion equation, 
the Gaussian solution of which is F s (qo,t) = exp(— q^Dt) = exp(— t /x a ). 

Simulation results for the mean square displacement of cations and anions at different temper- 
atures are displayed in Fig.|7J Subdiffusion in the intermediate time scale, another common feature 
of supercooled liquids, is quite pronounced, especially at low T . At long times, ion translations 
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Figure 7: Mean square displacement of (a) cations and (b) anions in the CGM description. 

tend to become normal diffusion. We notice that the time scale associated with the transition from 
non-Fickian to Fickian dynamics generally coincides with the structural relaxation time. 45 Another 
noteworthy feature is that for both A(f) and F s (qo, t), motions of PF 6 tend to be slower than those 
of EMI + . This is attributed to the fact that anions are heavier than cations. 
For additional insight, we have analyzed the relation between D and x a via 



where £ = 1 corresponds to the SE relation. We found that £ =0.87 and 0.92 for cations and 
anions, respectively (Fig. Ufa)). Thus the product Dx a develops a positive deviation from a constant 
value as T decreases (Fig. Ob)). This reveals a weak violation of the SE relation for our model 
RTIL system. 

The breakdown of the SE relation indicates that translational dynamics of the ionic liquid can 
not be described by the conventional diffusion equation, which is a continuity equation combined 
with a constitutive relation given by Fick's law. Specifically, it is assumed that the diffusion current 
is proportional to the spatial gradient of the particle concentration and obtains stationarity instan- 
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Figure 8: Breakdown of the SE relation in the coarse-grained ionic liquid, (a) Dashed and dotted 
lines are the fits with the scaling relation D ~ x a , where the exponent q is given by 0.87 and 0.92 
for cations and anions, respectively, (b) Dx a deviates from the constant value as the temperature 
decreases. Lines are drawn merely as a guide for the eyes. 

taneously in response to the external perturbation. Accordingly, the diffusion equation is valid 
only in the limit where the time is sufficiently coarse-grained to ensure instantaneous establish- 
ment of stationarity. If there is a significant delay in time before the system reaches stationarity, 
a crossover from the non-Fickian regime, characterized by anomalous diffusion and /3 relaxation, 
to the Fickian regime happens. We attribute the main cause of the delay to cage dynamics, i.e., 
immobile cages that last for a long time (see below). In this sense, the violation of the SE relation 
is a manifestation of dynamic heterogeneity, 32 which we turn to next. 



4 Correlated local excitations 

4.1 Decoupling of exchange and persistence times 

In ionic liquids, due to electrostatic interactions, a central ion is surrounded mainly by counte- 
rions in its immediate neighborhood, termed first solvation shell. Thus dynamics of the central 
ion will be influenced by those of the first solvation shell (cage) and vice versa. (In normal elec- 
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Figure 9: Displacement of a cation from simulation trajectories at (a) 300 K and (b) 600 K. Oc- 
currences of excitations at which the cation moves beyond d(=3.0 A) are marked with triangles. 
For a given initial position, the waiting time until the occurrence of the first excitation defines the 
persistence time T p , while the time interval between subsequent excitations yields the exchange 
times T x . 

trolytes, for example, this counter-ion cage, referred to as an ion atmosphere, tends to reduce the 
diffusion constant of the central ion at low concentration.) Suppose that the central ion undergoes 
thermal motions in the cage. Once in a while a large fluctuation enables it to escape the cage and 
subsequently the cage reorganizes. As the temperature is lowered, thermal fluctuations become 
suppressed and as a result, the frequency of the ion escape from the cage diminishes. While this 
picture may need quantitative elaboration, it is nonetheless useful to obtain a qualitative under- 
standing of "excitations" introduced below and their properties in glassy environments. 

We first introduce an excitation as a local event that an ion moves over a distance exceeding 
a threshold value d. The persistence and exchange times associated with excitations are then 
defined as follows: 22131 the persistence time T p is the waiting time t\ for an ion i to undergo its first 
excitation such that |r ; (?i) — r ; (0)| > d, and the exchange time T x includes a set of time intervals 
^2,^3, • ■ ■ between subsequent excitation events, i.e., \ri(t\ +t<i) —Ti(h) \ > d, |r,-(/i +^2 + ^3) —^i(h + 
h) \ > d, etc. Accordingly, the frequency of excitations gauges the mobility of ions. In the present 
study, we employ 3.0 A for d. We will return to discuss other possibilities for d later on. 
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Figure 10: Decoupling of persistence and exchange times for (a) cations and (b) anions in a coarse- 
grained ionic liquid. The persistence time T p and the exchange time T x are defined as the waiting 
times for an ion to produce a displacement exceeding d, taken as 3. OA, for the first time and 
thereafter, respectively. The solid and dashed lines correspond to the probability distributions of 
log T p and log T x , respectively. 



In a glassy environment where cage reorganization is slow, the likelihood of an ion undergoing 
a second excitation after its first one is higher than that of the initial excitation because the liquid 
structure disturbed by the first excitation generally provides a favorable environment for another 
excitation. In other words, the excitations are not governed by a Poisson process. In Fig.|9l typical 
trajectories of an ion at 300 K and 600 K are displayed. On each trajectory, the events of local 
excitations are marked with triangular symbols. As expected, excitations at 300 K are rare events; 
their occurrences are irregular and intermittent. This, for instance, leads to a non-exponential 
tail on the long-time end of the distribution of exchange times, i.e., time intervals between two 
consecutive excitations (see below). By contrast, the trajectory at 600 K is characterized by more 
regular and frequent excitations than that at 300 K. 

Figure [TOl exhibits the probability distributions of the logarithms of persistence and exchange 
times for cations and anions. At 800 K, the persistence and exchange times show nearly identical 



17 




1/800 1/400 1/300 
1/T (K" 1 ) 

Figure 11: Average persistence and exchange times versus the inverse temperature for d = 3.0 A. 
Lines are drawn as a guide for the eyes. 

distributions. This clearly indicates that the excitation events follow the Poisson process, viz., they 
are not correlated. As T decreases, two distributions become distinct and their difference increases. 
The center of the distribution for t p shifts to longer time much more rapidly than the corresponding 
distribution for T x . As a result, the deviation between the average exchange time (t x ) and average 
persistence time (t p ) rises markedly with lowering T (see Fig. [TT]). This decoupling of T x and 
T p observed here shows that the excitations become increasingly more correlated as T decreases, 
exposing the dynamically heterogeneous nature of our RTIL system at low T. We also notice 
dramatic enhancement of (t p ) at low T, compared to (t x ). For instance, (t p ) is longer than (t x ) by 
one order of magnitude at 300 K. This suggests the development of persisting immobile regions. 
According to a recent MD study in a similar RTIL, 46 the life time of stable contact ion pairs that 
seldom move a large distance as a pair can exceed a few nanoseconds. Such long-lived pairs could 
potentially be a candidate for immobile regions associated with long T p . Finally, as mentioned 
above, we notice that the distributions of both T x and r p develop a long non-exponential tail on 
the long-time end as T diminishes (note that the logarithmic time scale is employed in Fig. ITOl) . 
This is another indicator that the excitations at low T do not obey Poisson statistics and thus are 
correlated. 



18 



CO 



io 4 
io 3 

10 2 
10 1 
10° 



3Q0K + 

S50K 
400 K q 
475 K -m- 

6Q0K 
800K 



(a) 









- • 




,e' 
m' 










X) 


tzrV 




M' 



















□ 



CO 



1 



2 3 

d (A) 



4 5 6 




4 5 6 



Figure 12: Log-log plot of (a) the average exchange time (t x ) and (b) the average persistence time 
(t p ) for cations versus the threshold distance d. Lines are drawn as a guide for the eyes in (a) and 
(b). 

4.2 Analysis of threshold distance dependence 

In the previous subsection, we chose d = 3. OA as the threshold distance in the definition of exci- 
tations. The decoupling of the exchange and persistence times is present in a glassy environment 
with a physically meaningful value of d chosen, because it is attributed to strong correlations of 
local events. Here we examine how a different choice for the d value would influence our analysis. 
For example, if we choose a d value considerably less than the distance between neighboring ions, 
which is approximately 5 A for our system (Fig. [2]), excitations will correspond to small fluctua- 
tions of a central ion inside its counterion cage. We can easily imagine that such excitations seldom 
induce a considerable structural change in the local environment. Excitations exceeding ~ 5. OA 
on the other hand describe delocalized hopping or gradual drift, which is usually accompanied by 
irreversible structural changes. Thus, too large or too small a value for d in the working definition 
of excitations would not be able to capture, e.g., the decoupling of persistence and the exchange 
times even though the decoupling occurs irrespective of our choice of d. 

In Fig. [121 we display (t x ) and (t p ) of the cation versus d at various temperatures. We notice 
that (t x ) and (t p ) increase with d according to the power law, (t XiP ) ~ d s , where 8 ranges approx- 
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Figure 13: (a) Average exchange time (t x ) vs. average persistence time (t p ); (b) (t x ) vs. inverse 
diffusion constant D~ l ; (c) (t p ) vs. structural relaxation time T« for d = 5. OA, 3. OA, and v^A. 
Lines represent results fitted to (a) Eq. [6l (b) and (c) the scaling behaviors ~ x r , where x corre- 
sponds to D 1 and x a , respectively. In all cases, the value of the scaling exponents v and y, given 
for each d, is shown to increase with d. 

imately from 2 at 800 K to 4 at 300 K. In the scaling relation, the value 8 ~ 2 at 800 K corresponds 
to the diffusive regime, whereas larger values at lower temperatures indicate anomalous diffusion. 
Also the J-dependence of (t x ) and (t p ) becomes stronger as T decreases. 

In Sec. 13.21 the violation of the SE relation has been demonstrated via Eq.[5]with ^ =0.87 for 
cations. Likewise, the scaling relation 



(Tx) ~ (T P ) V , 



(6) 



is characterized by the exponent v, which is also smaller than unity. As shown in Fig.[T3{a), with 
d = 5.0 A, v is found to be 0.80 for cations. With d = 3.0 A, the corresponding v value is 0.67. We 
point out that the decoupling between (t x ) and (t p ) appears greater than that of the SE violation, 
that is, v is smaller than t, in Eq. [51 It is known that % a determined from F s (qo,t) is identified 
with the average persistence time, while the exchange processes contribute to the diffusion if d is 
comparable to the size of the ions (d > 5.0 A).— In this perspective, we analyzed the scaling 
behaviors of (t x ) with D 1 1 and (r p ) with x a . The results in Fig. fTST b) and (c) show that the former 
is characterized by weaker variations than the latter. We also notice that the coupling behavior 
becomes stronger with d, i.e., the exponent in the scaling relation increases, for both cases. But 
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Figure 14: Displacements of a cation and associated numbers of local excitation and brachiation 
events, m{t) and b(t). m{t) and b(t) were determined by counting the excitations and brachiation 
events for the central cation occurring during a time window of 500 ps along its MD trajectory. If 
neighboring anions revert back to original configurations within 200ps of their initial changes, the 



corresponding events were not considered as brachiation according to Ref. 148 



regardless of d, (t x ) shows a significant sublinear behavior in D 1 . This is due to the fragile 
characteristic of our model system, where the exchange events are correlated. 22 Therefore we 
expect that its scaling exponent would not reach unity even if a significantly larger value for d is 
employed in the definition of (t x ). By contrast, proportionality between (t p ) and % a is expected in 
the limit d — > 2n/qo. Note that the largest value of d employed in Fig. [T3f c) is 5. OA, which is still 
less than 2n/qo(= 7.32 A). 



4.3 Comparison between excitation and brachiation 

In an effort to understand ion displacements and related diffusion in RTILs, a mechanism based 
on brachiation processes, viz., a central ions moves by forming and breaking links to neighboring 
counterions through Coulomb interactions, has been proposed recently. 48 While this has some 
similarity to local excitations and cage dynamics, it is essentially a structure-based description 
that lacks dynamic correlation effects. To see this, we have performed a comparative analysis of 
brachiation and local excitations. In view of the brachiation length scale, — we reduced the cutoff 
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distance d slightly to v5A. In Fig. [H the MD results for a typical cation at 350 K together with 
the numbers of local excitations and brachiation events, m[t) and b(t), are displayed. We notice 
that while m(t) and bit) show a significant overlap in certain part of the trajectory, there is in 
general no strong correlation between the two. For instance, periods like 11ns < 13 ns and 
30ns < t < 36 ns are characterized by large b(t) but vanishing m{t). This means that neighboring 
ions reorganize in the presence of an immobile central ion. The opposite situation where m(t) is 
larger than b(t), e.g., 23ns < t < 25 ns and 27 ns <t < 29 ns, also occurs frequently, indicating ion 
translations without any considerable change of its neighbors. The weak correlation between local 
excitations and brachiation seems to suggest that the Coulomb interaction with neighboring ions 
does not play a major role in translational dynamics of individual ions in RTILs (see below). 

5 Role of Coulomb interaction 

Finally, we consider roles played by Coulomb interactions in the glassy behaviors of RTILs. 49 To 
this end, we compare with a model supercooled liquid that has dynamic characteristics similar to 
our EMI + PFg system but does not have long-range Coulombic interactions. 

We employ a model liquid system with the WCA potentials of Ref . 3 1 as a reference to quantify 
the influence of Coulomb interactions. For easy comparison, we follow Ref. [3l| to measure time 

1 /2 

and temperature in units of (m/O^/g/;) and £a/k^, respectively. We use the anion parameters, 
i.e., i = PF 6 , so that (m ( a?/e ; ,) 1/2 = 27.3 ps and £ U /k B = 180 K. Thus x a = 5030ps at T = 350K 
corresponds to X a = 184 at T = 1.94 in the new scaled units. The reader is reminded that the 
RTIL system at this temperature is characterized by strong nonexponential relaxation (Fig. [5]). By 
contrast, the WCA mixture does not exhibit glassy dynamics at all at this temperature. In fact, 
dynamics comparable to our RTIL at t = 1.94 occur at a much lower temperature t = 0.4 for the 
WCA system. 31 This reveals that at a given (scaled) temperature, structural relaxation dynamics 
in EMI + PF 6 in the CGM description are considerably slower than those in supercooled liquids 
characterized by short-range interactions only. This indicates that the Coulomb interactions in 
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Figure 15: Correlation coefficient px.y versus the inverse temperature, where X and Y denotes the 
instantaneous powers by the Coulomb force and the Lennard- Jones force. As the temperature is 
lowered, the negativity of linear correlations between X and Y becomes apparent. 

effect suppress structural fluctuations and enhance trapping in cage structures. 

To gain insight at the microscopic level, we briefly analyze instantaneous power inputs from 
the RTIL environment to an ion / via the Coulomb and LJ forces, i.e., Fp oul • V; and • v ; (v; = ion 
velocity). For simplicity, we consider only the united atom Tl for cations because it is the heaviest 
atom, located close to the cation center of mass. We examined the static correlation between the 
two instantaneous powers using the Pearson correlation coefficient px,Y'-— 

((X - (X))(Y - (Y))) 

PX.Y = , (/) 

<J X OY 

where X and Y denote Fp oul • v,- and F^ ■ v/, respectively, and Ox and ay are their standard devi- 
ations. The results in Fig. \T5\ disclose that power inputs arising from the LJ and Coulomb forces 
are anti-correlated. Therefore, LJ and Coulomb interactions play antagonistic roles in energy re- 
laxation of individual ions. If we calculate the time integration of Fp oul • v, and F^ ■ v,, however, 
the contribution from the former nearly vanishes and work to individual ions is delivered primarily 
by the LJ force(data not shown). This result together with our comparative analysis above paints 
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the picture that while the liquid structure of RTILs and thus the energy scale, i.e., f, relevant for 
glassy dynamics are mainly determined by Coulomb interactions, their relaxation is generally gov- 
erned by the LJ interactions. We believe this explains why dynamic aspects of glassy RTILs are 
very similar to those of non-ionic supercooled liquids despite their major difference in long-range 
interactions. We note that the anti-correlation of power inputs from the LJ and Coulomb forces and 
the dominance of the former in relaxation dynamics of the ionic liquid observed here also apply to 
vibrational energy relaxation in RTILs.— 

6 Conclusions 

We have introduced the coarse-grained model of ionic liquids and probed its dynamical behav- 
iors. Coarse-graining has simplified the geometry of the system and made dynamics accelerated, 
compared with the atomistic model. Nevertheless, the overall liquid structure and glassy dynamic 
properties such as nonexponential structural relaxation and subdiffusive behavior are preserved. 
Owing to the reduced number of atoms by coarse-graining, we have been able to perform extensive 
MD simulations at long times over a wide range of temperature, and investigated the temperature 
dependence of structural relaxation and diffusion. 

We found that our model for ionic liquids belongs to fragile glass formers, where x a exhibits 
strong non-Arrhenius dependence on the temperature. In addition, the SE relation is violated, 
which implies that diffusion is enhanced when compared with structural relaxation. We pay atten- 
tion to the apparent universality of the abovementioned dynamic features observed in a variety of 
glassy liquid systems, regardless of the nature of their molecular interactions. In previous studies 
of supercooled liquids, kinetic constraints have been successfully employed to explain the peculiar 
dynamic properties of supercooled liquids in view of facilitated dynamics. Dynamic facilitation 
emphasizes the dominant role of dynamic correlations in glassy dynamics rather than static prop- 
erties such as the structure factor and potential energy landscape. 51 In a similar fashion, we have 
defined a local excitation in the dynamics of our model and observe dynamic correlations between 
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them. Decoupling of persistence and exchange times has been shown to be highly correlated with 
local excitations. Such decoupling behavior exists regardless of the threshold distance d, which 
defines local excitations. However, d determines the physical meaning of local excitation events, 
especially at low temperatures. 

To understand the influence of the Coulomb interactions, we compared structural relaxation 
dynamics of the RTIL with those of non-ionic models of supercooled liquids. We found that glassy 
dynamics occur at a considerably higher temperature (in scaled dimensionless units) for the former 
than for the latter. We have also investigated instantaneous powers arising from the Coulomb and 
LJ forces. It was found that they are anti-correlated and their time integration is dominated by the 
latter. These results seem to indicate that relaxation dynamics of RTILs are dominated by the LJ 
interactions, while the Coulomb interactions exert a strong influence on the liquid structure and 
thus set the temperature scale for glassy dynamics. 

At low temperatures, immobile regions persist for a long time due to sparse excitations. Once a 
local excitation occurs, subsequent displacements of ions are more probable, which tends to intro- 
duce mobile regions. Thus, decoupling of persistence and exchange times is a plausible explanation 
for dynamic heterogeneity of glassy liquids. In the future, we plan to investigate spatiotemporal 
correlations of local excitation events in order to characterize in more detail dynamic heterogeneity 
in the RTIL. 
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